# Note!! You will need to run 02 - main models.R in the global environment
# Otherwise, the code won't run, since the bootstraps also use the old models and data

# Load packages and install them if necessary (automatic)
need = c(
  "dplyr",           #data structuring/cleaning
  "janitor",         #data structuring/cleaning
  "sandwich",        #robust se's
  "scales",          #data structuring/cleaning
  "stringr",         #data structuring/cleaning
  "tidyr"            #data structuring/cleaning
)     
have = need %in% rownames(installed.packages())
if(any(!have)) {install.packages(need[!have])}
pack <- lapply(need, library, character.only = TRUE)
# Clean environment
rm(have, need, pack)

# Start with setting a seed for replication purposes (just in case, also set in function)
set.seed(28)

# Optional but handy (turn off scientific notations)
options(scipen = 999) 

# Load functions for bootstraps and computing output
source("D:/dv_rep/Source/boot_functions.R")
# Load data
source("D:/dv_rep/Source/dat_gss.R")

#### Bootstrapped models ####

# Set N bootstraps to 1000
n_bootstraps <- 1000

# Bootstraps without controls (note: together, these bootstraps will take a considerable amount of time, up to several hours!)
genDec_boot_gssPolviews  <- perform_bootstrap(gss_Polviews, n_bootstraps, gen = "genDec", include_controls = FALSE)
genDec_boot_gssSwelfare  <- perform_bootstrap(gss_Swelfare, n_bootstraps, gen = "genDec", include_controls = FALSE)
genDec_boot_gssRacgender <- perform_bootstrap(gss_Racgender, n_bootstraps, gen = "genDec", include_controls = FALSE)
genDec_boot_gssCultmoral <- perform_bootstrap(gss_Cultmoral, n_bootstraps, gen = "genDec", include_controls = FALSE)
genDec_boot_gssImm       <- perform_bootstrap(anes_Immig, n_bootstraps, gen = "genDec", include_controls = FALSE)

# Bootstraps with controls
genDec_cboot_gssPolviews  <- perform_bootstrap(gss_Polviews, n_bootstraps, gen = "genDec", include_controls = TRUE)
genDec_cboot_gssSwelfare  <- perform_bootstrap(gss_Swelfare, n_bootstraps, gen = "genDec", include_controls = TRUE)
genDec_cboot_gssRacgender <- perform_bootstrap(gss_Racgender, n_bootstraps, gen = "genDec", include_controls = TRUE)
genDec_cboot_gssCultmoral <- perform_bootstrap(gss_Cultmoral, n_bootstraps, gen = "genDec", include_controls = TRUE)
genDec_cboot_gssImm       <- perform_bootstrap(anes_Immig, n_bootstraps, gen = "genDec", include_controls = TRUE)

# Extract parameters with extract_short function
param_genDec_gssPolviews  <- extract_short(genDec_boot_gssPolviews, gss_Polviews)
param_genDec_gssSwelfare  <- extract_short(genDec_boot_gssSwelfare, gss_Swelfare)
param_genDec_gssRacgender <- extract_short(genDec_boot_gssRacgender, gss_Racgender)
param_genDec_gssCultmoral <- extract_short(genDec_boot_gssCultmoral, gss_Cultmoral)
param_genDec_anesImm      <- extract_short(genDec_boot_gssImm, anes_Immig)

cparam_genDec_gssPolviews  <- extract_short(genDec_cboot_gssPolviews, gss_Polviews)
cparam_genDec_gssSwelfare  <- extract_short(genDec_cboot_gssSwelfare, gss_Swelfare)
cparam_genDec_gssRacgender <- extract_short(genDec_cboot_gssRacgender, gss_Racgender)
cparam_genDec_gssCultmoral <- extract_short(genDec_cboot_gssCultmoral, gss_Cultmoral)
cparam_genDec_anesImm      <- extract_short(genDec_cboot_gssImm, anes_Immig)

# Make data frames for each model
tab_genDec_gssPolviews     <- make_Tab(param_genDec_gssPolviews)
tab_genDec_gssPolviews$var <- "Ideological self-placement"
tab_genDec_gssSwelfare     <- make_Tab(param_genDec_gssSwelfare)
tab_genDec_gssSwelfare$var <- "Social welfare"
tab_genDec_gssRgender      <- make_Tab(param_genDec_gssRacgender)
tab_genDec_gssRgender$var  <- "Race and gender"
tab_genDec_gssCmoral       <- make_Tab(param_genDec_gssCultmoral)
tab_genDec_gssCmoral$var   <- "Culture and morality"
tab_genDec_anesImm         <- make_Tab(param_genDec_anesImm)
tab_genDec_anesImm$var     <- "Immigration"

# Models with control variables
ctab_genDec_gssPolviews     <- make_Tab(cparam_genDec_gssPolviews)
ctab_genDec_gssPolviews$var <- "Ideological self-placement"
ctab_genDec_gssSwelfare     <- make_Tab(cparam_genDec_gssSwelfare)
ctab_genDec_gssSwelfare$var <- "Social welfare"
ctab_genDec_gssRgender      <- make_Tab(cparam_genDec_gssRacgender)
ctab_genDec_gssRgender$var  <- "Race and gender"
ctab_genDec_gssCmoral       <- make_Tab(cparam_genDec_gssCultmoral)
ctab_genDec_gssCmoral$var   <- "Culture and morality"
ctab_genDec_anesImm         <- make_Tab(cparam_genDec_anesImm)
ctab_genDec_anesImm$var     <- "Immigration"

# Combine to a single data frame
tab_genDec_gssDat <- plyr::rbind.fill(
  tab_genDec_gssPolviews,
  tab_genDec_gssSwelfare,
  tab_genDec_gssRgender,
  tab_genDec_gssCmoral,
  tab_genDec_anesImm
) 

# Add column to specify that these models are without control variables
tab_genDec_gssDat$controls <- paste("Without")

# Models with controls
ctab_genDec_gssDat <- plyr::rbind.fill(
  ctab_genDec_gssPolviews,
  ctab_genDec_gssSwelfare,
  ctab_genDec_gssRgender,
  ctab_genDec_gssCmoral,
  ctab_genDec_anesImm
) 

# Repeat
ctab_genDec_gssDat$controls <- paste("With")

# Combine to single data frame
tab_genDec_gss <- plyr::rbind.fill(
  tab_genDec_gssDat,
  ctab_genDec_gssDat
)

# Save output
save(tab_genDec_gss, file = "D:/dv_rep/Datasets/genDec_gss.Rdata")

# Make tables
tab_genDec_gss <- clean_Tab(tab_genDec_gss)

tab_genDec_gss <- tab_genDec_gss %>% 
  mutate(
    p_star = case_when(
      p < 0.001 ~ "***",
      p < 0.01  ~ "**",
      p < 0.05  ~ "*",
      TRUE ~ ""
    )) %>%
  mutate(across(c(coef, se), ~ as.numeric(.))) %>%
  mutate(across(c(r2, r2_dev, adj_r2), ~ as.numeric(format_fixed_decimals(., decimals = 3)))) %>%
  mutate(coef_p  = paste0(format_fixed_decimals(coef, 2), " (", format_fixed_decimals(se, 2), ")", p_star)
  )

# Convert to wide format
tab_Wide_genDec <- tab_genDec_gss %>%
  select(term, controls, coef_p, facet, boot, var) %>%
  pivot_wider(., names_from = var, values_from = coef_p) %>%
  janitor::clean_names() %>%
  select(
    term,
    facet,
    boot,
    controls,
    ideological_self_placement,
    social_welfare,
    race_and_gender, 
    culture_and_morality,
    immigration,
  ) %>%
  as.data.frame()

rows_Tab <- tab_genDec_gss %>%
  select(controls, facet, boot, var, r2, r2_dev, adj_r2, N) %>%
  distinct %>%
  mutate(N = as.numeric(as.character(N))) %>%
  pivot_longer(cols = c(r2, r2_dev, adj_r2, N), names_to = "term", values_to = "coef_p") %>%
  mutate(coef_p = ifelse(term %in% c("r2", "adj_2", "r2_dev"), format_fixed_decimals(coef_p, 3), coef_p)) %>%
  mutate(coef_p = as.character(coef_p)) %>%
  mutate(coef_p = case_when(term == "r2" & boot == "Boot" ~ paste0(coef_p, " (", lead(coef_p), ")"), TRUE ~ coef_p)) %>%
  filter(term != "r2_dev") %>%
  pivot_wider(., names_from = var, values_from = coef_p) %>%
  janitor::clean_names()

# Create a list of unique groups in tab_Wide
unique_groups <- unique(tab_Wide_genDec[c("controls", "facet", "boot")])

# Initialize an empty data frame to store the combined result
combined_df <- data.frame()

# Iterate over unique groups and append rows_Tab after each group in tab_Wide
for (i in 1:nrow(unique_groups)) {
  group <- unique_groups[i, ]
  group_rows_Tab <- rows_Tab %>%
    filter(controls == group$controls, facet == group$facet, boot == group$boot)
  
  group_tab_Wide <- tab_Wide_genDec %>%
    filter(controls == group$controls, facet == group$facet, boot == group$boot)
  
  combined_group <- bind_rows(group_tab_Wide, group_rows_Tab)
  combined_df <- bind_rows(combined_df, combined_group)
}

# Mean table
meanTab_genDec <- combined_df %>%
  filter(facet == "Mean" & boot == "Boot") %>%
  select(!boot & !facet)

meanTab_genDec <- meanTab_genDec %>% 
  pivot_wider(., names_from = controls, values_from = c(
    "ideological_self_placement",
    "social_welfare",
    "race_and_gender",
    "culture_and_morality",
    "immigration"
  )) %>%
  mutate_all(~ replace_na(., "")) 

# Write to csv
write.csv(meanTab_genDec, file = "meanTab_genDec.csv", row.names = FALSE)

# Divergence
divTab_genDec <- combined_df %>%
  filter(facet == "Divergence" & boot == "Boot") %>%
  select(!boot & !facet)

divTab_genDec <- divTab_genDec %>% 
  pivot_wider(., names_from = controls, values_from = c(
    "ideological_self_placement",
    "social_welfare",
    "race_and_gender",
    "culture_and_morality",
    "immigration"
  )) %>%
  mutate_all(~ replace_na(., "")) 

# Write to csv
write.csv(divTab_genDec, file = "divTab_genDec.csv", row.names = FALSE)

# Sorting
sortTab_genDec <- combined_df %>%
  filter(facet == "Sorting" & boot == "Boot") %>%
  select(!boot & !facet)

sortTab_genDec <- sortTab_genDec %>% 
  pivot_wider(., names_from = controls, values_from = c(
    "ideological_self_placement",
    "social_welfare",
    "race_and_gender",
    "culture_and_morality",
    "immigration"
  )) %>%
  mutate_all(~ replace_na(., "")) 

# Write table to csv
write.csv(sortTab_genDec, file = "sortTab_genDec.csv", row.names = FALSE)
